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Метод побудови математичних моделей 
складних процесів на засадах 
генетичних алгоритмів 


Розроблений метод побудови математичних моделей складних процесів на засадах генетичних 
алгоритмів. У порівнянні з багаторядним алгоритмом групового урахування аргументів даний метод 
дає можливість синтезувати математичні моделі будь-якої складності без попереднього вибору числа 
рядів селекції. Метод знайде застосування для підвищення точності прогнозів коливних процесів, 
наприклад, рівня ріки Дністер. 


Вступ 


Цілому ряду процесів економічного та екологічного характеру притаманна гармо- 
нічна складова з некратними частотами, яка моделюється наступним виразом: 


С(г)- Аз з У; зітіоо ; В; созісо .)), (1) 
узі 


де г - такти відліку часу, / - І, 2, 3,..5 М; А, 4;, В; - параметри гармонічного ряду 
(5 о; сфрудо; - некратні частоти, / - 11,2, 3,.... 

Для того, щоб можна було б оцінити параметри ряду (1), необхідно виконання умо- 
ви (1) У 2 Зт- І. Ідентифікація параметрів моделі (1) відбувається у декілька етапів |21. 
Перший етап - обчислення балансових коефіцієнтів о.,, із умови мінімізації нев'язки 
М-т 


во У, 2) 


їзтя 
-І 
де Б, -х(ічт)- У а,(я(і я р) ві - р) ві -т), іст 1,М -т; є(0) - відліки реаліза- 
ре0 


ції процесу у моменти часу, що симетрично розміщені відносно довільної точки 7. 
Отже, будемо розв'язувати задачу 


з М-т т- І 1 
тіп : До)з 2 Зіт 7 З ві з (3) 
- ізтч ре0 

де б-(ад Ор), - вектор вагових коефіцієнтів; вій зді т) Щі -т); Ре здічрічняі -р); 


Т - символ транспонування матриць. 
Задачу (3) запишемо у матрично-векторній формі 


пий й Ла)- (2, -ка) -Биа) , (4) 
єОт-- 1) вИ) Жтч 1) Дт-- 2) Я (т) М зт) 202) 
дрова «От 2)--8(2) р Жт- 2) Дт. 3) -Я(т-- Ї) за «От 1) 90) 
вх) в(м «т мні АУ-тч 1) нім-т-) у м-н дет 
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Мінімізація виразу (4) приводить до нормального рівняння Гауса, яке у матричній 
формі матиме такий вигляд: 


ТрообТз 
ро Сем ій Ду Р (5) 
Із останнього рівняння можна знайти 
І 
ре т Па 
а-(кТЕ, ЕЇРи: (6) 


Використовувати формулу (6) можна лише тоді, коли розмірність вектора а 
невелика і матриця Е/ Е, є добре обумовленою. Якщо така умова не виконується, то 


для знаходження с слід розв'язувати рівняння (5) одним із числових методів, 
наприклад, методом Гауса зі зворотним ходом (31. 
т-1 
На другому етапі складається рівняння ад У а, со5 (ро) - соз (ту) для 
раї 
довільної частоти с», яке за допомогою рекурентного співвідношення (1 со8( ро) - 
- 2 со5(( р -1) е)сов а - со5((р - 2) о), р з 2,т приводиться до алгебраїчного рівняння 
т -то степеня відносно со5 аз 
Ре" Раз Рач Р «0, (7) 
де г асо. 

Суть третього етапу у тому, що, знаючи вагові коефіцієнти а.,, можемо скласти 
рівняння (7), розв'язок якого відносно 2 дає змогу однозначно визначити частоти гар- 
монік о, / - І,т. Тепер задача полягає в оптимальному синтезі гармонічного ряду (1). 
Відомі два підходи |2| до вирішення поставленої задачі. Перший з них передбачає 
викреслювання гармонік у різних комбінаціях із повного ряду ;, / 2 іт. Другий 
метод грунтується на ідеях багаторядних алгоритмів групового урахування аргументів 


(МГУА). Суть якого в наступному. Деяка реалізація вихідної величини процесу дов- 
жиною М розбивається на три частини: М, - навчальну послідовність, М, - переві- 


рочну 1 М є - екзаменаційну. На першому ряді селекції за всіма заданими точками виді- 


ляються всі можливі тренди гармонічного ряду; максимальне число трендів (1| тудх 2-- 


Із них вибирається не єдиний тренд, а 4 трендів, які у найбільшому степені задо- 
вольняють вибраному критерію селекції. Після цього обчислюється 4 залишків (залиш- 
ком називають різницю ординат коливного процесу 1 кожного із трендів першого ряду). 
На другому ряді селекції із кожного залишку знову виділяється 4, трендів. Із всієї мно- 
жини 44, отриманих трендів другого ряду за тим же критерієм селекції вибирається 4 
кращих трендів цього ряду 1 т.д. Число рядів селекції збільшується до тих пір, поки 
зменшується величина критерію селекції. На останньому ряді селекції вибирається єди- 
ний розв'язок, який відповідає мінімуму критерію селекції. 

Недоліком першого підходу до вирішення поставленої задачі є необхідність пере- 


то. 
бору великого числа варіантів, яке визначається як сума 5, - У, С, - 2" - 1, що потребує 
іні 


значних затрат машинного часу. Для другого підходу характерним є те, що у результаті 
реалізації багаторядного алгоритму МГУА неможливо отримати математичну модель у 
явному вигляді. 
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Синтез математичної моделі гармонічного ряду 


Запропонований інший підхід побудови математичних моделей коливних про- 
цесів, який базується на ідеях генетичних алгоритмів. Вся реалізація ово величини 
процесу або явища розбивається на три частини у такій пропорції |2| М, -0,7М, 
М от 0,2М 1 М, - 0,1М. Для множини даних М, - М о ВИЗНАаЧАаЮТьЬСЯ вагові коафійнах 
ти а, З рівняння (4) за методом виключення Гауса з вибором головного елементу |З. 
Розв'язок рівняння (7) відносно змінної 2 дає можливість знайти частоти с;, / «І т. 
Тоді на множині точок М,, -- М, необхідно знайти параметри моделі (1) 40, 4,, В, |4). 
Тому утворимо упорядковану структуру (хромосому) довжиною м, в якій на і-му місці 
буде стояти нуль або одиниця в залежності від того, чи частота с»; вилучена із ви- 
браного повного ряду т, чи залишена. Набір хромосом утворює популяцію. У задачі 


синтезу моделей коливних процесів функцією пристосованості, яка дозволяє вибрати 
найбільш пристосовані особі з популяції, виступає комбінований критерій селекції 51 


релі є 87, (8) 


де пі - критерій зміщення, який обчислюється за такою формулою: 
М 
2 
У (є (Ю)- є, (5) 
пі - З! У , (9) 
Хі 
із! 
де В - функція нев'язки, що визначається як (2); 2. (К), є; (5) - величини, значення 


яких обчислені відповідно на множині точок М за формулою (1), а коефіцієнти моделі 
(1) знайдені відповідно на множинах М, Мо і Му. 


Генетичний алгоритм складається з наступних кроків ||. 

КІ. Формування початкової популяції (ініціалізація). На першому кроці роботи 
алгоритму випадковим чином формується популяція з І осіб, кожна із яких є хромо- 
сомою, довжиною т. 

К2. Оцінка пристосованості хромосоми у популяції. Для кожної хромосоми об- 
числюється критерій селекції (8). У відповідності з моделлю (1) формується матриця 


1 5іп Ф, со5 Ф, 5іп Ф, со5 Ф, ж 5 Ф,, со5 Ф, 
й 1 віп(29,) | со8(22,)  зіп(29,)  сов(20,) ... віп (2), ) со (20, ) 
1 зіп(Мо,) со8(Мо,) зіп(Ме,) сов(Мо,) .. зіпП(Ма,) сов (Мо, ) 


У хромосомі подвоюємо одиниці і нулі тому, що кожній частоті и Й відповідає 
пара коефіцієнтів 4;, В;. Оскільки у моделі (1) завжди присутній коефіцієнт 40, то 
до хромосоми на першу позицію добавляємо одиничний ген. У відповідності до 
сформованої хромосоми СЛИ;д із матриці Є формуємо нову матрицю ЕК", шляхом 


пеу 


вилучення тих стовпців із матриці Б, які асоційовані із нулями хромосоми СП). 
Матрицю Ек, утворюють перші М, - М о стовпці матриці Е,,,, а другу - останні Мя. 
На множинах точок М, 4 Мо і МУ; обчислюються ненульові коефіцієнти 4,, 4, В, 
моделі (1) шляхом розв'язку нормального рівняння Гауса 
Ек Екк Ре» 10) 
Е5 Ру Аз з Бо 85, а1) 
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де Ак «Я 0, ія) вія) А вію, 2 449), 405), 809), 45) ві), ді8) ві? г 


ті" ті 


вектори параметрів моделі, яка асоційована з черговою хромосомою із початкової по- 


дою), яз -(80.0.805) - 


вектори експериментальних даних на множині точок М; -н Мо і Му. 


пуляції, обчислені за формулами (10) 1 (11); як Чай. 2. 


За відомою сукупністю коефіцієнтів А, і А; моделі (1) на множині точок М 
обчислюють «(в)- розро зо «(5)- Кона Значення критерію селекції знаходять для 
кожної хромосоми за формулою (9) і в результаті отримують множину значень 0), 


і-1,М ,де М - кількість хромосом у популяції. 
КЗ. Перевірка умови зупинки алгоритму. Визначають 
Ртіп - тіпр). (12) 
М 


Є 
Якщо мінімальне значення (12) критерію селекції (9) не перевершує деякого 
додатного значення є, то відбувається зупинка алгоритму. Зупинка алгоритму також 
може відбутися у випадку, коли його виконання не приводить до покращення функції 
пристосування, або у тому випадку, коли алгоритмом уже виконано задане число 


ітерацій. Після виконання однієї із трьох умов із популяції вибирається хромосома СА", 
для якої ВИКОНУКЕРНЯ умова (12). Після операції ПОДВОЄННЯ і приєднання одиничного 


гена отримуємо -- Сйто . Ця хромосома задає структуру моделі оптимальної складності і 


формує матрицю ХЕ" таким чином, що із початкової матриці ЕК вилучаються стовпці, 
які асоційовані з нульовими генами хромосоми Со . Перерахунок параметрів моделі 
(1) здійснюється на множині всіх точок початкового масиву даних. 

КА. Селекція хромосом. Найбільш поширеними методами селекції |6| є метод 
рулетки і метод турнірної селекції. У даному алгоритмі використано турнірний метод. 
При турнірній селекції всі хромосоми розбиваються на підгрупи, з наступним вибором 
із кожної утвореної підгрупи хромосоми з найкращою пристосованістю. 

К5. Формування нової популяції потомків здійснюється за допомогою двох ос- 
новних операторів: схрещування і мутації. Оператор мутації з ймовірністю Р, змінює 
значення гена в хромосомі на протилежне, тобто з І на 0 чи з 0 на 1. Оператор схре- 
щування складається з двох етапів. На першому етапі вибирається найкраща хромосома 
з підгрупи осіб за критерієм селекції. У результаті отримуємо нову популяцію хро- 
мосом, до якої застосовують оператор другого етапу - оператор схрещування. 

Із пулу родичів М(Кк) випадковим чином з ймовірністю Р. утворюють пари. 
Генерується випадкове число Р. із інтервалу (0; 1) і якщо його значення не більше, ніж 
Р. то над парою хромосом здійснюється схрещування, інакше пара хромосом 
залишається без зміни. Потім для кожної пари родичів розігрується позиція гена (локус) 
в хромосомі, яка визначає точку схрещування /,,. Фіксація точки схрещування зводи- 
ться до випадкового вибору цілого числа із інтервалу |/; І, - Ї. Дія оператора схре- 


щування приводить до того, що із пари родичів утворюється нова пара потомків. 
Далі відбувається перехід до К2. 


Приклад синтезу математичної моделі 
з використанням генетичного алгоритму 


Часто паводки призводять до розмивання берегів, що може викликати пошко- 
дження газопроводу «Торжок-Долина» (Твано-Франківська область), який прокла- 
дено біля русла р. Дністер. Тому для забезпечення надійної та безаварійної роботи 
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газопроводів, що проходять через ріку Дністер, важливе значення має прогнозування 
його паводків. За період з 01.04.07 по 31.08.07 велись спостереження за рівнем води 
р. Дністер у районі с. Нижнів Івано-Франківської області. 

Аналіз зміни рівня води у р. Дністер за вказаний період засвідчує, що з часом 
має місце тренд /(г), який має лінійний характер, та існує гармонічна складова Сг), 
зумовлена сезонною зміною метеорологічних умов |71, тобто 

НН, Сі) НИ), (13) 
де Н, - поточний рівень води, см; С(г) - гармонійна складова рівня води, см; М(г) - 
лінійний тренд, см. 

Із залежності (13), яка визначається зміною рівня води у р. Дністер, був виді- 
лений лінійний тренд 

кг) є во виг, 4) 
де 00, 0; - параметри лінійного тренду, які знайдені за методом найменших квадратів. 


Результат роботи програми, яка написана з використанням розробленого методу у 
середовищі Маїї аб, відтворює рис. 1, на якому знаком «0» відмічені експериментальні 
дані, а «Р» - результат розрахунку за формулою (1). Після виділення із експерименталь- 


них даних гармонічного тренду, отримали залишок - Н, - Н рн С(о)чн пі). Величина 
Н, є функцією параметрів, що визначають погодні умови у районі спостережень, тобто 

Н; ЗТ, Лов Ла Лов Р) (15) 
де Т, - середньодобова температура повітря, С; /, - кількість опадів, мм/добу; м, - сере- 
дньодобова швидкість вітру, м/с; р, - середньодобовий барометричний тиск, мм рт. ст., 
г - поточний дискретний час, К - 1, 2,3 - зсув у часі. 


40 - т т ; т т - - й 


2015-4188 


Рівень, см 


Рисунок І - Гармонічний тренд коливного процесу (р. Дністер) 


Будемо розглядати ділянку ріки як деяку систему, що характеризується сукуп- 
ністю вхідних величин Хх - (х),Х»,зоху) і вихідною величиною у (рис. 2). У нашому 


випадку - х; - Т,, Хо, З |» Хз З Дт» Ха З Дощ Хо хо Ух Я р УЗН. 


Система у 


Рисунок 2 - Структурна схема системи «ділянка ріки - спостерігач» 


Співвідношення (15) будемо шукати у вигляді полінома 


М-І К 5 її 
А- ХУ а Пх/ з (16) 


іс | узі 
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де М -- кількість членів полінома; а; -- коефіцієнти полінома; 5 ;; -- степені аргументів, 


п 
які повинні задовольняти обмеженню - 5,5); Хт. Число членів М полінома (16) ви- 


-І 
значають за такою формулою (21: 4 


- (тп), рег (7) 

Вхідні величини х), /- 1,Е, які є аргументами виходу системи у, у кожному спо- 
стереженні / приймають певне значення так, що їх сукупність утворює матрицю Х. 
Допустимо, що нам відомі параметри а;, і- 1,М - І моделі (16). Тоді за відомими 


значеннями величин хі) можна обчислити 


МАК убує 1 
уза Уа Пп зле (18) 
із | узі 
Систему рівнянь (18) зручно подати у матрично-векторній формі 
у-Ка, (19) 


де уУз(уруУзноекМ у - обчислене значення виходу моделі (16) у кожній точці 

спостережень; ЕК - матриця розміром М х М , елементи якої добутки аргументів при 
ре т Й б 

параметрах а;; а - (ад,а,. ам!) - вектор параметрів моделі (16). 


Знаючи 7, і у,, і - 1,М., можна обчислити критерій апроксимації 
М 
ее 2 
ла Мо: (20) 
І-І 
Мінімізація критерію (20) і приводить до співвідношення 
ЕТКа- ТУ, (21) 


яке називають нормальним рівнянням методу найменших квадратів (МНК). 
Безпосередньо з рівняння (21) можна знайти 
а «|ЕТЕ'ЕТТ. 02) 

Використовувати формулу (22) можна лише тоді, коли розмірність вектора 
параметрів а невелика і матриця Е"Е є добре обумовленою |8). Якщо така умова не 
виконується, то для розв'язку рівняння (22) слід використовувати один із числових 
методів, наприклад, метод Гауса з вибором головного елементу |З. 

На вихід системи у накладається перешкода е, так що спостерігачу доступна 
тільки величина У - у -е,. Якщо допустити, що е адитивна і має нормальний закон 
розподілу, то оцінки параметрів моделі (18) є незміщеними та ефективними (38). 

Як правило, структура моделі (18) невідома, що приводить до необхідності до- 
вільного вибору як числа функцій, так і вигляду самих функцій у моделі (18). Тому для 
вибору структури моделі (18) був запропонований індуктивний метод самоорганізації 
моделей |1), ідейну сторону якого визначає теорема Геделя. Стосовно задачі визначення 
структури моделі (18) геделівський підхід означає застосування зовнішнього критерію, 
який дає можливість однозначного вибору єдиної моделі із заданого класу моделей. От- 
же, всі дані, які отримані у результаті експерименту, розбиваються на дві частини М, - 
навчальна і М, - перевірна. 


Для вибору структури моделі використовують критерії регулярності 
Мо 
2 (У РА у 
Л2 (0 ) - 41 хо 
УТ; 
ізі 


(23) 
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і мінімуму зміщення 
М 


у (1)-у, (0) 
д7(в,0)- 1! (24) 


М 
УЗ 2 
ї«І 


Якщо вибраний критерій регулярності (23), то вибирають наступний розподіл да- 
них експерименту (2) М, -0,7М і М, -0,3У, а при виборі критерію (24)- М, - 0,5М 


і М 20,5М. 


Реалізація індуктивного методу самоорганізації моделей здійснюється поетапно: 
перший етап - генерація моделей-претендентів (у певному порядку підвищення склад- 
ності); другий етап - відбір найкращої моделі за критерієм селекції (23) або (24). 

Розрізняють три способи генерації моделей-претендентів. Перший із них - ком- 
бінаторний метод, який вибирає найкращу модель на основі одного із критеріїв се- 
лекції із сукупності моделей, що отримують із виразу (13) шляхом прирівнювання до 
нуля деяких його коефіцієнтів. Другий спосіб відомий як метод групового врахуван- 
ня аргументів (МГУА), у якому генерація моделей здійснюється на основі багаторядної 
процедури. Третій метод подібний до другого, але на кожному ряді селекції часткові 
моделі утворюють шляхом прирівнювання до нуля певного числа їх коефіцієнтів. 

Недоліком комбінаторного методу є необхідність перебору великого числа моде- 
лей, що практично неможливо реалізувати за допомогою сучасних ЕОМ. МГВА поро- 
джує моделі, у яких фігурують проміжні змінні кожного із рядів селекції, що значно 
утруднює процес переходу до вхідних змінних системи, що моделюється. Це відносить- 
ся і до третього методу. 

Найкращим є комбінаторний метод, оскільки він дає можливість отримати 
модель, де аргументами виступають вхідні величини системи. Для зняття проблеми 
великої розмірності застосуємо генетичний підхід. Як емпіричну модель розглянемо 
поліном (18) степені т, де частина параметрів приймає значення нуль. Утворимо 
упорядковану структуру довжиною М , в якій на і -му місці буде стояти одиниця або 


нуль в залежності від того, чи параметр а,, і-0,М - І моделі (18) відмінний від нуля, 
чи нульовий. Отже, із початкової популяції хромосом шляхом еволюційного відбору 
потрібно вибрати таку хромосому, яка забезпечує найкраще значення функції присто- 
сування (мінімальне значення критерію селекції (23) або (24)). Алгоритм розв'язку 
поставленої задачі аналогічний раніше розробленому для виділення гармонічного 
тренду з тією лише різницею, що не здійснюється операція подвоєння генів у хро- 
мосомах. 
Знайдені залежності /(г), С(/) і у, дають можливість знайти 


в-СО)нВИ)чу,, (25) 


де С(г), (г) і у, - обчислювались відповідно за формулами (1), (14) 1 (19). 


Графік залежності (25) показаний на рис. 3, де «-» позначені обчислені значення 
за формулою (25), а значком «о» відмічені експериментальні значення рівня води у 
р. Дністер. Із графіка бачимо, що мають місце досить задовільні збіги між розрахун- 
ковими і експериментальними даними. Адекватність моделі перевірялась за допомогою 
коефіцієнта кореляції Ку, - 0,98746 , що свідчить про високий степінь кореляції між 


величинами У; 1 у). 
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Метод побудови математичних моделей складних процесів... 2г 


Рівень, м 
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Рисунок 3 - Залежність рівня води у р. Дністер від погодних умов 


Висновки 


Таким чином, застосування ідей генетичних алгоритмів до побудови матема- 
тичних моделей дає можливість вибрати оптимальну за структурою адекватну модель і 
значно зменшити об'єм обчислень. Остання обставина відкриває широкі можливості для 
побудови складних моделей як фізичних явищ, так і технологічних процесів. Ефектив- 
ність алгоритму підтверджена на прикладі моделювання зміни рівня води у р. Дністер у 
залежності від погодних умов. Отримана модель може бути використана при прогно- 
зуванні повеней, що є досить актуальною задачею. 


т равура 


1. Ивахненко А.Г. Индуктивньй метод самоорганизации моделей сложньх систем: |монография| / 
А.Г. Ивахненко. - К. : Наукова думка, 1981. -- 296 с. 

2. Ивахненко А.Г. Справочник по типовьтм программам моделирования / А.Г. Ивахненко, Ю.В. Коппа, 
В.С. Степашко и др. ; под ред. А.Г. Ивахненко - К. : Техніка, 1980. - 180 с. 

3. Вержбицкий В.М. Основьт численньх методов : |учебник для вузов| / В.М. Вержбицкий. - М. : 
Вьсшая школа, 2002. - 840 с. 

4. Горбійчук М.І. Індуктивний метод побудови математичних моделей газоперекачувальних агрегатів 
природного газу / М.І. Горбійчук, М.Ї. Когутяк, Я.І. Заячук // Нафтова і газова промисловість. - 
2008. -- Хо 5.- С. 32-35. 

5. Ивахненко А.Г. Помехоустойчивость моделирования : |монография| / А.Г. Ивахненко, В.С. Степа- 
шко. - Киев : Наук. думка, 1985.-216с. 

6. Рутковская Д. Нейронньке сети, генетическиє алгоритмьт и нечеткиєе системьі / Д. Рутковская, М. Пи- 
линьский, Л. Рутковский ; |пер. с польск. Й.Д. Рудинскогої. - М. : Горячая линия- Телеком, 2004. -- 452 с. 

7. Горбійчук М.І. Математичні моделі прогнозування стоку р. Дністер для запобігання техногенних 
аварій магістральних газопроводів / М.І. Горбійчук, О.В. Пендерецький // Науковий вісник Івано- 
Франківського національного технічного університету нафти і газу. - 2008. -- Ме 2 (18). - С. 30-33. 

8. Ермаков С.М. Математическая теория оптимального зксперимента / С.М. Ермаков, А.А. Жигляв- 
ский : (учеб. пособиг). - М. : Наука, 1987. - 320 с. 


М.И. Горбийчук, М.А. Шуфнарович 

Метод построения математических моделей сложньх процессов на основе генетических алгоритмов 
Разработан метод построения математических моделей сложньжх процессов на основе генетических 
алгоритмов. По сравнению с многорядньм алгоритмом учета аргументов данньй метод даєт 
возможность синтезировать математические модели любой сложности без предварительного вьтбора 
числа рядов селекции. Метод найдет применение для повьшения точности прогнозов колеблющихся 
процессов, например, уровня реки Днестр. 
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